What are the elementary excitations of the BCS model in the canonical ensemble? 
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We have found the elementary excitations of the exactly solvable BCS model for a fixed number 
of particles. These turn out to have a peculiar dispersion relation, some of them with no counterpart 
in the Bogoliubov picture, and unusual counting properties related to an old conjecture made by 
Gaudin. We give an algorithm to count the number of excitations for each excited state and a 
graphical interpretation in terms of paths and Young diagrams. For large systems the excitations 
are described by an effective Gaudin model, which accounts for the finite size corrections to BCS. 
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I. INTRODUCTION 

The paradigmatic model to study the superconducting 
properties of metalsi and nuclei^ is the pairing model pro- 
posed by Bardeen, Cooper and Schrieffer. The ground 
state (GS) and excitations of the BCS model are well 
known in the grand canonical (g.c.) ensemble, and ex- 
plain the behavior of systems with large number of par- 
ticles. However for small systems, such as nuclei or 
nanograins, one is forced to work with a fixed number 
of particles, where the g.c. BCS wave function, includ- 
ing its projected version, are not adecuate. The problem 
is due to the strong pairing fluctuations, which a mean 
field approach cannot deal with properly. An alterna- 
tive approach is provided by exact numerical methods, as 
the DMRG^, but their complexity somehow obscures the 
Physics behind. Fortunately enough the reduced BCS 
model, characterized by a unique pairing coupling g, is 
exactly solvable, as was shown by Richardson*. This ex- 
act solution has been recently used in connection with 
superconducting nanograms^ (see Ref. for a review) . 

Most of the exact studies deal with the GS and the ex- 
cited states that are obtained by breaking Cooper pairs. 
However one must also consider the promotion of pairs 
to higher energy levels (bosonic pair-hole excitations). 
This paper focus on the later type of excitations since 
the former ones can be easily included into our formal- 
ism. We shall indicate the peculiar dynamics and the 
unusual counting properties exhibited by the excitations 
of the exactly solvable BCS model, some of them with 
no analogue in the standard picture of Bogoliubov quasi- 
particles. These features account for the exact finite size 
corrections to the thermodynamic limit, obtained from 
the standard BCS treatment. 

In Section [H] we introduce the superconducting sys- 
tem described by Richardson and Gaudin models, and we 
present here our conjecture about what represents their 
elementary exciations. Section ITTT1 is devoted to classify 
the excited states according to what will be interpreted 
as elementary excitations in Section llVI In this respect 
a diagramatic representation of the excited states turns 
out to be very useful. The thermodynamic limit of our 
theory is presented in Section lTvl which confirms our con- 



jecture about the elementary excitations. A comparison 
with the BCS theory is also presented. Finally, SectionlVl 
contains our conclusions. 



II. RICHARDSON AND GAUDIN MODELS 

Let us consider a fermionic system with N single parti- 
cle energy levels. The reduced BCS Hamiltonian decou- 
ples the levels which are singly occupied and one is left 
with those that are either empty, |0) (hole) or doubly oc- 
cupied, &]|0) (pair) with an energy Ej (b^ = cj, + c] _ is a 
hard core boson operator that creates a pair of two time- 
reversal states). We shall suppose that the singly occu- 
pied levels have been removed. Since the latter decou- 
ple, their effect can be considered easily by adding their 
free energy Sj to the total energy. The complete system 
will be treated elsewhere. The reduced BCS Hamiltonian 
reads, 

JV TV 
3=1 Jij'=l 

where G is a dimensionful coupling constant. The stan- 
dard model employed to study nanograins is given by 
the choice £j — d(2j — N — 1), where d = oj/N is the 
single particle energy level spacing and ui/2 is the De- 
bye energy^. The coupling G can be written as G = gd, 
where g is dimensionless. The value of the bulk BCS 
gap, AbcSi °f the equally space model is given by 
Abcs — A/2, where A = aj/sinh(l/g). 

The eigenstates of (Q with M pairs are given by4 




where the parameters {E fl }^ l L 1 satisfy the Richardson 
equations 




,u=l,...,M. (3) 
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FIG. 1: Real part of for the equally space model with 
M = 7V/2 = 20 pairs and N G = 0, 1, 2 excitations. 

The total energy of the state © reads E = X^u=i 
The number of solutions of Richardson Eqs. © is given 
by the binomial coefficient = ( M ), and coincides 
with the dimension of the Hilbert space, H M , of states 
with M pairs distributed into N different levels. Then it 
is natural to label Richardson states by a set of integers 
/ = {ji , . . • , 3m } corresponding to the value that the pair 
energies E^ take at g = 0, i.e. some of the e,-'s. 

The BCS model can be mapped into a spin system 
which at g — > oo has SU(2) symmetry. Based on this 
fact Gaudini made the conjecture that given a solution 
of Eqs. and taking the limit g — » oo, 

a subset of them, say {E a (oo)}a=i, remain finite and 
satisfy the equations 

N 1 N G 

° = E— -f-E^-f ' ° = 1 -->Afc. (4) 

while the remaining M — A*g roots tend to infinity and 
satisfy Eqs. 10 with all £,-'s set to zero. The number Nq 
of finite roots takes values from to M. The number of 
solutions of Eqs. (0J is given by g?at g = C^ G — C^ G _ 1 7 . 
Therefore in the large g limit the Richardson's so- 
lutions would be classified in terms of Nq according 
to Table ||| Consistency is guaranteed by the equation 

U M — 2^N G =0 aN o- 



# of solutions 


d M 


d,M-i 
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do = 1 


Efj, finite 
Ey,, infinite 


M 



M - 1 
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1 
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TABLE I: Classification of roots in the g — > oo limit. 

We show in this paper that Gaudin finite energies rep- 
resent the elementary excitations of the superconducting 
system in the canonical ensemble. Their peculiar disper- 
sion relations and the unusual counting properties will 



FIG. 2: Position of the M = 20 pairs of the states of fig. □ 
at g = 1.5. The arcs Ti 1 (19 pairs) and Ti 2 (18 pairs) are a 
slight modification of the GS arc Fi (20 pairs). 



account for the finite size corrections to the mean-field 
BCS treatement of superconductivity. 

This result is motivated by the excitation energy for 
large values of g, namely, E exc = E — Eqs ~ guiNcW — 
(Nq — 1)/N], and the gap A ~ gui. Thus, in the large N 
limit the excitation energy goes as E exc ~ NqA. Which 
allows us to think of the state as a set of Nq elementary 
excitations contributing each with an energy A to the 
total energy. 

In Section Hvl we extend this result for the whole range 
of g. In the meanwhile next Section is devoted to the 
proof of Gaudin's conjecture given in Table U We also 
obtain a formula which gives the number of finite Gaudin 
energies Nq for a given Richardson configuration J, and 
therefore the number of elementary excitations. 



III. CLASSIFICATION OF EXCITED STATES 
A. Simple examples 

Let us first consider the simplest examples given by 
the excited states with one and two energies remaining 
finite, i.e. No — 1 and 2. Representatives of these, to- 
gether with the GS are shown in Fig. ^ which depicts to 
the real part of the energies, and in Fig. |3 which shows 
the distribution of the energies in the complex plane for 
g = 1.5 and a system with M = 20 pairs at half filling, i.e. 
N = 2M. As a general feature we see that for small g all 
parameters E^ are real, and as g grows some of them col- 
lapse and become complex conjugate pairs, which share 
their real part (this corresponds to two curves merging 
into a single one in Fig. Fig. [21 shows how the ener- 
gies E^ arrange themselves into an arc which opens up 
to infinity as g — > oo. 

The state of Fig. corresponds to the GS of the sys- 
tem, and it is labelled by Iq = {1, 2, . . . , M}, which at 
g = is identical to the Fermi state (FS). As g — > oo 
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all the roots become complex and scape to infinity. Ac- 
cording to Table nj this is the only state where this may 
happen, hence N G (I ) = 0. 

The lowest excited state I\ = {1, . . . , M — 1, M + 1} is 
shown in Fig. The last root Em, which is equal to 
Sm+i at g = 0, stays finite as g — *■ oo, while the remaining 
M - 1 roots go to infinity thus N G {h) = 1- All the 
states with N G = 1 can be obtained from the FS by: 
1) promoting the nearest pair below the Fermi level (FL) 
into one of the N — M empty levels above it, or 2) moving 
the nearest hole above the FL into one of the M occupied 
levels below it. The state 1\ can be obtained in both 
ways. Hence the number of N G = 1 excited states is 
di = (N — M) + M - 1. 

The state of Fig. [TJi, J 2 = {1, . . . , M - 2, M, M + 2}, 
has Ng(I2) — 2. All the states with N G = 2 can be ob- 
tained in three different ways from the FS: 1) moving the 
two pairs just below the FL into the N — M empty levels 
(C^~ M states), 2) moving the two holes just above the 
FL into the M occupied levels (C^ 1 states), or 3) moving 
one of the M — 1 pairs in the FS, except the closest to 
the FL, into one of the N — M — 1 vacancies above the 
FL, except the closest to the FL ((N - M - 1)(M - f) 
states). The state with two holes just below and two 
pairs just above the FL is generated by the rules 1 
and 2, thus the number of states is the expected one, 
d 2 = Cz~ M + Ci 1 + (N — M — l)(Af - 1) - 1. This 
example shows that the value of N G for a generic state 
depends dramatically on the arrangement of holes and 
pairs around the FL. 

B. N G (I) formula 

We now turn to the evaluation of Nc(I) for a general 
state. One naively expects that this formula should be 
given by the sum of pairs, N p , and holes, Nh, above and 
below the FL respectively, i.e. Ng(I) = N p + Nh- In fact 
N p = Nh since every pair above the FL comes from a hole 
below it. However this ansatz does not always work as 
we have already seen above. For example, according to 
this formula, the state I\ of Fig.^J) would have 2 instead 
of Nq = 1, while the state I2 of fig-H: has N G = 2, which 
is the correct value. 

Let us introduce for convenience the occupation rep- 
resentation of the states /, where a pair, a hole and the 
FL are depicted as •, o and | respectively. In the cases 
discussed above we obtain: 7 = »- -- »»»|ooo---o, 
Ii = »- -- »»o|»oo---o, and /2 = , ''' ,0, | 0,0 '''°- 

We have found an algorithm to compute Nq {I) ■ Given 
an integer £ > 0, let us split / into three disjoint sets, I — 
A( U Be U Ci, where Ag contains the lowest M — £ levels, 
Bi the next 2£ levels and Ce the remaining N — M — £ 
ones. For I = 0, the set Bq is empty, and Aq (resp. Co) 
contains all the levels below (resp. above) the FL, while 
for £ > 1 the set Bi contains the nearest £ levels above 
and below the FL. As an example, let us choose a state 

of the form ^3 = • - V - ■ •oo»«»«|o»«o ■ • • o. For £ — 2 the 



partition of ^3 is given by {• ■ ■ ■ •oo»»}{»»|o»}{»o ■ • ■ o}. 

Let us define the number of pairs and holes for each set, 
i.e. N* h (X = A e ,B e ,C e ). The algorithm giving N G (I) 
is 

N G (I) = min* = o 1 ... l2 jv p N G {I,£) (5) 
N G (I, £) = + min(JVjf < , ) + N^ . 

Applying this formula to Is one gets {N G (I 3 , £)}f— — 
{4,5,4,3,4} and thus, N G (h) = 3. The value of £ u 
which minimizes N G (I,£), is given in this case by 3 (in 
general £j is not equal to N G (I)). The result of this 
formula is bounded, N p < N G < 2N p , and therefore 
any state with a finite N p would contain a finite number 
of Gaudin energies. The uncorrelated counting formula 
proposed earlier coincides with the case £j = 0, since 
N G (I, 0) = N p + N h (notice that N h = Np and N p = 

The physical mechanism underlying Eq. JSJ is the col- 
lective behavior of the holes and pairs that occupy the £j 
closest levels to the FL. In a certain sense, lj measures 
the range of correlation involved in the creation of the el- 
ementary excitations out of an initial pair-hole configura- 
tion. However, this correlation can be lifted introducing 
a shifted Fermi level FL(£i) defined by moving the FL an 
amount of £1 levels downwards (resp. upwards) whenever 
N^ 1 is lower or equal (resp. greater) than Np . This new 
Fermi level defines a new Fermi state out of which the 
excited state with N G finite energies is obtained by the 
creation of uncorrelated pairs above and holes below the 
FL(£j). This construction provides a pathway to the g.c. 
formulation as discussed later. 

The formula presented here allows us to prove Gaudin's 
conjecture by looking for all the states with a given N G , 
and finding out that its number corresponds to djv G as 
stated in TableQJ the same way we already did with those 
with N G = 1 and 2. 

C. Young diagrams representation 

The correlated behavior of the excitations is made 
more explicit by a pictorial representation of the states. 
The idea is to associate to every set / a path 7/ with 
N links on the square lattice Z 2 , starting at the ori- 
gin (0,0). This is achieved by associating a horizontal 
link directed to the right, to every hole o, and a ver- 
tical link directed upwards, to every pair •. The map 
starts from the lowest energy level and ends at the high- 
est one. For example the path associated to the Fermi 

M N—M . , ...... . 

state lo = ••••• I o • • ■ o is a polygonal line joining the 
points (0, 0) -»• (0, M) —* (N — M, M). If / describes a 
state with M pairs and N energy levels, then the path 
7/ ends at the point (A'' — M, M). The number of these 
sort of paths is C M , which is precisely the dimension of 
the Hilbert space Wj^- 

In Fig. [SJi we depict the occupation and path represen- 
tations of the state ^3 which yields N G = 3, in agreement 
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FIG. 3: a) The path and Young diagram of I3. b) Real part 
of Ey, for 73. For g large enough there is a real root (1) and a 
complex root (2). 



with the numerical results shown in Fig. [3J3. Moreover, 
Fig. illustrates the fact that any state I gives rise to 
a Young diagram (YD) Yr, whose boundary is formed 
by the links which belong either to 7/ or to 7/ , but 
not to both. The YD of the Fermi state is by construc- 
tion empty, i.e. Y/ = 0. These YD's capture the basic 
properties of the excitations. First of all, Ng(I), given 
by Eq. J5J, coincides with the number of squares on the 
longest SW-NE diagonal on Yj (see Fig. Et). This fact 
provides a geometrical meaning to Na(I) and leads to 
a combinatorial proof of Gaudin's conjecture, which can 
be stated as follows: d^ G is the number of YD's, Yj, as- 
sociated to the paths 7/ , which have Nq squares on their 
longest SW-NE diagonal. The proof of this conjecture 
uses the methods of Ref. |^. This result serves to classify 
the excitations in terms of YD's. For example the states 
with Nq = 1 and 2 discussed above correspond to the 
YD's shown in Fig. 0] 



IV. THERMODYNAMIC LIMIT 

As we explained in the previous Section, by increas- 
ing g, M — Nq of the energies become complex and 
arrange themselves into an arc which escape to infinity 
for large g, while the remaining Ng stay finite with their 
positions barely modified (see Figs.[IJ),c andOl. 

Following the procedure presented in Refs. 
we take the large N limit keeping M/N, g and Ng finite. 
In this limit the arc formed by the energies in the com- 
plex plain becomes dense, and allows for a continuous 
formulation. In particular, the GS corresponds to an arc 
Ti in the complex energy plane, which in the g — » 00 
limit goes to infinity. 

Excited states contain finite energies in addition to the 
arc. A given finite root E a can be either real or complex. 
In the former case we shall call it a 1-string. In the latter 
case E* a is also a root, which together with E a form a 2- 
string (an example of such states can be seen in Fig. [3b). 
There are also 3-strings formed by one real root and two 
complex ones, having approximately the same real part, 
and so on. In general {E a }^°i is a combination of strings 
with several lengths. The remaining M — Ng roots con- 
dense into an arc Tj, which is a slight perturbation of the 
GS arc IV . In Fig. [5] we depict T] and Tj for the two 
excited states I\ and I2 shown in Figs. nj),c. 

Taking into account these considerations, and using 
the methods of Refs. PIUH one can show in the large N 
limit that the excitation energy of a Richardson state / 
is given by 



N G 



E exc = J2 V(E a - £o) 2 + A 2 



(G) 



a=l 



where e$ is twice the chemical potential, and the energies 
E a satisfy the modified Gaudin equations 



JY 



1 



E 



x R( Sj )( Sj - E a ) fe a R(E P )(E - E a ) 



(7) 




FIG. 4: YD's corresponding to Ng = 1 and 2. 



Another properties of these diagrams are: i) the pair- 
hole transformation of the states induces a transposition 
of their associated YD's, ii) the main NW-SE diagonal 
on a YD coincides with the FL (see Fig.|3^), and ii) the 
number of boxes of Yj is the excitation energy of / (in 
units of 2d) at g — for the equally space model. 



with R(E) = ^(E-e ) 2 + A 2 . As g — > 00 one has 
A ~ gu and Eqs. Q become Eqs. |@J|. 

The excitation energy given by Eq. © fits quite well 
the excitation energies of our prototype example (N — 
40, M = 20), as shown in Fig.El This also exhibits the 
linear behavior of the excitation energy for g — > 00, i.e. 
Eexc ~ Ng&- as stated in Section^ and in full agreement 
with the large g behavior of Eq. @ . Thus we can extend 
our conjecture to the whole range of g. Namely, any ex- 
cited state is composed by Nq elementary excitations as- 
sociated to the finite Gaudin energies. The Hilbert space 
spanned by these excitations has therefore a dimension 
d,N G = ( Na ) — (jv^-i)- Hence, it is reasonable to call this 
new type of excitations Gaudin pairs or gaudinos. 

In order to compare our results with the BCS stan- 
dard solution lets consider the excitation energy given by 
a real Cooper pair in the Bogoliubov approach, which is 
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FIG. 5: Excitation energies E axc = E — Egs < 14d for 
M = 20 pairs at half filling. There are 44 =13+26+5 states 
corresponding to Ng = 1,2,3 respectively. The particle-hole 
symmetry reduces these numbers to 25 = 7+15+3. 



given by Je\ + A 2 (notice that A = 2A B cs), and span 



a Hilbert space of dimension 



The standard Bo- 



goliubov quasiparticle with an energy ye 2 + A 2 /2 would 

have to be compared with excitations involving broken 
Cooper pairs. Since E a in Eq. © lies between two en- 
ergy levels, with £j+i — £j = 2d ~ 1/N (e.g. in Fig. 2Jj 



E 20 (oo) = with e 20 < E 20 < £21), E a 



and d 



N G 



I N 



0(l/N), 



N-Ng+1 

consistent within 0(1 /N) corrections, as it is well known 
from the existing relation between a canonical and a 
grand canonical ensemble formulation in the statistical 
physics. 

It is important to notice how the Bogoliubov excita- 
tions are uncorrelated with respect to the BCS ground 
state. We already pointed out that the correlation 
present in our formulation is lifted by choosing a shifted 
Fermi level FL(£j). This Fermi level is within a dis- 
tance 0(1 /N) from the original one. The selection of 
a new Fermi level leads a new Fermi sea (with a differ- 
ent number of particles), allowing for a grand canonical 
formulation in a natural way. 

In summary, our gaudinos will yield the same results 
as the BCS theory in the extrict N — 00 limit, and will 
account for the exact corrections to the bulk results for 



Ng 



—]. Therefore, our theory is 



the finite size superconducting grains for all the physical 
observables and thermodynamic properties. 



V. CONCLUSIONS 

We have shown in this paper that the elementary exci- 
tations of the exactly solvable BCS model in the canoni- 
cal ensemble can be explained by the Gaudin model and 
have no counterpart in the Bogoliubov picture of quasi- 
particles. Their peculiar dispersion relation and the un- 
expected counting properties, which are due to the cor- 
related behavior of pairs and holes around the Fermi 
level, provide the exact finite size corrections to the BCS 
bulk results, valid for large systems in the g.c. ensem- 
ble. These excitations, together with those obtained by 
breaking Cooper pairs, supply the complete spectrum of 
the canonical BCS model. A formula to compute the 
number of elementary excitations for any given state was 
also proposed. 

We explained how the description in terms of gaudinos 
agrees with the Bogoliubov picture in the thermodynamic 
limit to leading order in N. In the case of broken pairs, 
which was not presented here, the mechanism is identical. 
It is of interest to study how the phase of the supercon- 
ducting order parameter emerges from this fixed number 
of particles formulation. It will be intimatelly related to 
the possibility of choosing a shifted Fermi level in the 
large N limit (which looses the correlation of the exci- 
tation), allowing the introduction of ground states with 
different number of pairs. 

Although we used as an example a system of equally 
spaced levels, the results are more general, and apply to 
any distribution of levels. This assertion is based on nu- 
merical calculations considering broken Cooper pairs. In 
this case bloked levels are removed, and we are left with 
a non-equally spaced spectrum, obtaining again the same 
general results. In the case of a non-constant pairing we 
also expect the qualitative picture presented here to hold. 
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